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ABSTRACT 



I'M 

Methods of investigating internal waves from bathythermograph 
information are discussed. A series of 651 hourly observations 
taken at Ocean Station "P" is analyzed by two methods to determine 
the characteristics of the internal waves present. One method 
involves a study of the fluctuations in depth of six selected 
isotherms; the other involves a study of the fluctuations in 
depth of the top, center, and bottom of the thermocline. An 
objective method is used to determine the depths of these 
features from the bathythermograph data. 

Spectra of these features are calculated by Fast Fourier 
transform, hanned, and smoothed, and these spectra are compared 
with previously published results obtained from a different 
method of analysis of the same time series. 
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INTRODUCTION 



Tidal period internal waves in the deep ocean were observed 
as long ago as 1907, by the Swedish oceanographer 
Pettersson 1909 . Little work was done in describing these 
waves or investigating their cause before 1950, although Stokes 
developed the theoretical treatment which could describe internal 
waves- -neglecting geostrophic effects-- in 1847. This theory 
was used by Ekman jl904 in explaining the ’’dead water" 
phenomenon in Norwegian fjords. Although internal waves have 
been observed to exist for more than half a century, only 
recently have mathematical tools and equipment been developed 
which enable investigators to more closely analyze these waves in 
terms of periodic components and energy distribution. Although 
oceanographers are now in a position to say a great deal about 
the various frequencies which are present in a spectrum of 
internal waves, the presence of many of these frequencies 
and the causes of them have not yet been explained. 

The most useful tools for analysis of internal wave data 
are Fourier analysis and power spectrum analysis. Fourier 
analysis, although a tool of mathematicians for decades, 
required great expenditure of time to be carried out before the 
advent of high-speed digital computers, and has only recently 
been brought into use in analyzing internal wave records. Power 
spectrum analysis, developed and refined by Blackman and 
Tukey ^1959^ has been used by electronic and communications 
engineers for several years, but has only recently been applied 
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to ocean wave analysis. Programming either of these methods for 
use on a digital computer makes possible analysis of large amounts 
of data accurately in a reasonably short time. 

In order to describe internal waves, it is at least necessary 
to obtain a time or space series of observations for analysis. 

The most common method of obtaining a time series is to observe 
temporal changes in thermal structure at a fixed location. This 
is sometimes done from weather observation ships located at deep 
ocean stations. One of these. Ocean Station Tt P n , located in the 
northeastern Pacific Ocean at 145 ° West Longitude, 50° North 
Latitude and maintained by the Department of Transport of the 
Dominion of Canada, has been making bathythermograph observations 
at a maximum interval of 12 hours since July 1952 . From time to 
time, bathythermograph observations are made from the ship 
occupying Station ”P” at four hour intervals, or more frequently. 
The resulting bathythermograph data provides an excellent time 
series for analysis of internal waves. 

The data for use in this thesis consisted of hourly 
bathythermograph observations, obtained with mechanical bathy- 
thermographs from 2000 GMT, 5 August 1961 to 2200 GMT, 

1 September 1961, a total of 651 hourly observations. The data 
were furnished in digital format on punched cards by the 
National Oceanographic Data Center, Washington, D. C. Digital 
analysis was performed on an IBM 360 computer at the Naval 
Postgraduate School Computer Facility, Monterey, California. 
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Ba thy thermograms of the data used in this thesis may be found 
in Manuscript Report Series No. 106, published by the Fisheries 
Research Board of Canada, Programmed by the Canadian Committee on 
Oceanography . 
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II. OBTAINING DATA FOR POWER SPECTRUM ANALYSIS 



One of the most commonly used methods of obtaining data for 
analysis in describing internal waves has been following the 
depth fluctuations of selected isotherms. This method was used 
by Tabata jl965j at Ocean Station ,T P" in the northeastern 
Pacific Ocean, and by Seiwell j*1965j in the North Atlantic Ocean. 
Although this method has been quite successful in the past, its 
validity depends on rather specific conditions. Internal waves 
propagate along isopycnal surfaces, and have their greatest 
amplitude at a density discontinuity. An isothermal surface 
will closely approximate an isopycnal surface for propagation 
of an internal wave only if the water column is nearly isohaline 
and the wave is of small amplitude . In this case one can approx- 
imate the change in density due to salinity and the change in 
density due to pressure fluctuations by zero. In choosing the 
isotherms of interest, Tabata selected isotherms which were 

within the permanent thermocline, assuming that the thermoeline 

( 

would oscillate as a unit, although phase differences between 
thermoclines would not be indicated in the power spectrum. 
Seiwell, in the North Atlantic, chose the 20° C isotherm as being 
close to defining the top of the thermocline, or depth of the 
mixed layer. 

Defant |l932j in examining the data collected by METEOR 
anchor stations, used a similar method, although he examined 
temperature fluctuations at fixed depths. In order to compare 
his results with that of other investigators, one must select 
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a depth near the top of the thermocline, since other investigators 
used isotherms near the top of the thermocline for analysis. 
Decreases in temperature would indicate the approach of a crest, 
and increases in temperature the approach of a trough. 

LaFond and Rao 195*+j obtain their raw data from the mixed 
layer depth, disregarding any changes above or below this depth, 
which corresponds to the top of the permanent thermocline. This 
would appear to be a more realistic approach, since the top of 
the thermocline will approximate a density discontinuity, and 
internal waves should have their greatest amplitudes at this 
point, and thus be most easily detected. 

In this thesis, the fluctuations in depth of the top, center, 
and bottom of the thermocline were used to investigate internal 
waves. The depths of these features were obtained using a method 
originated by Boston 1966i in which the thermocline is made to 
approximate a Gaussian distribution of depths and related 
temperatures. The depths thus obtained were subjected to a 
Fourier analysis, and the coefficients obtained were harmed 
and smoothed over discrete frequency bands. 

For comparison of results of isotherm analysis, the depths 
of six selected isotherms were followed through the time series 
and subjected to analytical procedures identical to those used 
in the spectral analysis of the thermocline. 



9 



III. OBTAINING THERMOCLINE DEPTHS 



A. THE GAUSSIAN THERMOCLINE 

A gaussian, or normal, distribution of temperature (T) as a 
function of elevation (Z) is given by the equation 



T(Z) = 



V2tF 



~ 2 “. 

e dz 



-(Jo 



where Z is positive upwards. The frequency function corresponding 
to this distribution function is 



dT-Z) 1 -4 

dZ V 2tt 



Using a basic statistical approach, defining a normalized variable 

7 , Z-m 

^ - a 

where m is the mean of the Z values, and CT is the standard 
deviation of the array of Z values, we obtain a new frequency 
function which is 



[ Z-nV)2 

dT(Z) _ 1 
dZ cjV2tt 

This function has the following characteristics: 

1. It is symmetric about the point Z = m 

2. It has two symmetric points of inflection at Z = m-CT 

3. It has a maximum rate of change at Z = m-CT~\f3 

A change in the numerical value of m causes a displacement of the 
curve in the vertical direction, but does not alter its form. 
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A change in has the effect of altering the scale in the 
vertical direction. The smaller the value of , the more 

concentrated the curve is about the .point Z = m. 

If the water surface is defined to be Z = 0, and the depth of 
the water considered to be infinite, then the function defining 
the vertical temperature distribution is 



t(Z')= 



a V2'rr 




(Z-mr 

*5 



dZ 



-oo 



The center of the curve T(Z') versus Z is found at the point 
Z = -m. 

If this curve is considered as a bathythermograph trace, the 
following definitions may be made. The center of the distribution 
of T(Z') versus Z is the center of the thermocline, called Z Q . 

The top and bottom of the thermocline will correspond to the 
points of maximum rate of change of slope of this distribution. 

If Z t is used to designate the depth of the top of the thermocline 
and Zj-jf- to designate the depth of the bottom of the thermocline, 
these definitions become 

Z c = -m 
= -m+cr 
z bt = -m-cr 

These features are illustrated in Figure 1A. 

To obtain the necessary values to apply these definitions, 
consider the first and second moments of the (T,Z) curve about 
the origin, Z = 0. The nth moment of the curve about the origin 
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n — 

is, by definition, the mean value of the array Z , called Z , 
where n is a positive integer. With a normal curve, the first 
moment is identical with the mean, which corresponds to the 
center of the thermocline. 

T(0) 

Zc' / z < T > dT 

T(-«) 

where T(0) is the temperature at the surface and T(-OO) is the 
temperature of the bottom of the lower isothermal layer which, 
so far, has been assumed to be the bottom of the ocean. The 
second moment 

T(G) 

^ = TiO-T(-°o)/ z(TfJT 

f(-«) 

can be used to compute O , since the variance of the distribution, 
'J - , is defined as 

cr^= Z~-m~=Z^ -Zq 

In practice, several approximations to this exact solution 
must be made. The temperature T(-JO) is measured at the bottom 

of the bathythermograph trace. The differentials, dT, become 

C- 

finite differences, ^T, and the integrals are replaced by 
summations . 

In changing from an integral representation to a summation 
over finite intervals, a certain amount of error in the value 
of is introduced, due to the fact that summation in this 

fashion is not a continuous process, whereas integration is. 
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To compensate for this loss of accuracy, a weighting factor must 
be introduced in all computations involving CT . Grosfils jl968 
calls this factor sk, and bases the size of the factor on the 

'i 

size of OT. 



For St^1.0 C sk = 1.30 

For St<0.5° sk = 1.47 

Intermediate values of St will have intermediate values of sk, 
but these values were chosen as being those most commonly 
encountered. Using this factor, the equations defining the top 
and bottom of the thermocline become 

Z-j. = -m + C7 sk 
-m CT sk 



Using this approach to define the thermocline, one can fully 
describe a thermocline knowing the temperature at the surface, 
the temperature at the bottom of the bathythermograph cast, and 
depths corresponding to some constant temperature interval. 

The center of the Gaussian thermocline can now be defined as 






and since 



then 



Z(Tn)i 



1 



T 1 ‘ T bt 



Z(T n ) 



cr 2 = Z(T n )f-Z 



2 

C 



Using these factors, the top and the bottom of the Gaussian 
thermocline may readily be determined. 
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In the event the distribution of depth with temperature is 
not Gaussian, the inaccuracies built into this objective 
definition are intolerable. Whether or not this distribution 
is Gaussian can be determined by computing the third moment of 
the (T,Z) curve about the origin. This, however, is unnecessary 
when one applies a non-Gaussian approach to the thermal structure. 

B. THE NON-GAUSSIAN THERMOCLINE 

The approach that was used by Boston and by Grosfils was 
based on depth data given at constant temperature intervals. The 
data used for this thesis were temperature data given at 
constant depth intervals. This difference does not affect 
the treatment of the data except for the determination of the 
temporary center of the thermocline. The details of the method 
used in this thesis to find the center of the thermocline are 
given in Section IV. 

When the depth of the center of the thermocline is found, by 
whatever method the data requires, the data is then processed as 
two Gaussian thermoclines — one based on the depth distribution of 
temperatures above and including the depth Z Q , the other based on 
the depth distribution of temperatures below and including the 
depth Z c . 

In order to locate Z^, only the upper portion of the (T,Z) 
curve is considered. This portion of the curve is "reflected” 
about Z c , to produce the curve shown in Figure IB. From this 
curve C~~ is computed by methods previously described and, 
using this value, Z-j- may be calculated. 
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c 

FIGURE I 
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In order to determine Z^, a similar procedure is followed, 
except that in this case, the lower portion of the (T,Z) curve 
is "reflected" about Z c , to produce the curve shown in Figure 1C. 
Note that it is possible for this resulting curve to include 
negative depths, but since the depths are squared in 

O 

computing CT , this does not result in a narrowing of the 
distribution through an underestimate of CT . 

It can be seen that in the event the thermocline is Gaussian, 
treating it in this fashion will not alter the results obtained 
by treating it as Gaussian. Thus, every (T,Z) input may be treated 
as a non-Gaussian distribution, regardless of its true form. 

This is the method which was used in processing the ba thy thermogram 
data for this thesis. 
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IV. DIGITAL METHODS USED TO OBTAIN DEPTH INPUTS 



The bathythermograph data was obtained from the National 
Oceanographic Data Center in digital format. For each 
bathythermograph trace, the digital data was contained on three 
or four computer data cards, including the sequential number of 
the observation, the date and hour of the observation, 
temperature in degrees Centigrade at each five meters of depth, 
and such other bits of synoptic information as the location of 
the observation, the ship making the observation, and a number 
identifying the cruise during which the data was taken (Figure 2) 

The first phase of data reduction was to transfer the 
sequential number of the observation, the time and date, and 
the temperatures at the first twenty depths--which correspond to 
temperatures from the surface to a depth of 95 meters--to a 
computer data card. This was not essential for data reduction, 
but was done for convenience as it reduced the size of the data 
deck by three quarters. By inspection of ba thy thermograms for 
the entire series of observations the permanent thermocline was 
in all cases determined to be shallower than 95 meters in depth. 

The data deck thus obtained was then analyzed to obtain the 
depths of the top, center, and bottom of the thermocline using 
a method based on a digital computer program developed by 
Grosfils | 1968 . Grosfils' program, however, was written for 

data consisting of depths corresponding to equally spaced 
temperatures, and since the NODC format differs from this, some 
modifications had to be made in the program. 
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FIGURE 2 

N0DC Data Format 



The computer first reads the data from the data deck, then 
calls subroutine INTERP, which fits straight lines through the 
given data points, then determines the depths of temperatures at 
intervals of 0.1°C from the surface to 95 meters. The temperature 
is read by the computer as an integer — in tenths of degrees — and 
left in this form for the interpolation process, to eliminate 
truncation error in manipulating the temperatures. 

In the program written by Grosfils, the center of the 
thermocline was defined as the point where the second difference 
of depths was greater than zero. Due to the format of the data 
used in his program, however, a great deal more information was 
known about the thermocline than is given in the NODC format. 

For this reason, the following method was used to define the 
center of the thermocline, utilizing the array of temperatures 
before interpolation. 

In subroutine INTERP, the first differences of temperatures 
were calculated, from the surface (corresponding to T^) downward. 
The ith first difference is defined as 

ATi = T. - T i+1 

The first time the ith first difference is less than -0.3°C, T^ is 
stored as ITTOP and the computer then begins calculating first 
differences again at a depth of 45 meters. When the ith first 
difference then becomes greater than -0.2° C, the computer stores 
Tj_ + -l as ITBOT, The temperature of the center of the thermocline 
is then defined as the arithmetic mean of ITTOP and ITBOT, and is 
stored as ITC. The computer then interpolates depths and 
temperatures as desc rived previously. 
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Once the interpolation is completed, subroutine PURGE is 
called to insure that only one depth is recorded for each 
temperature value, which would not be the case if any temperature 
appeared more than once in the data which was read in. PURGE 
examines the interpolated temperatures, and in the event that 
any temperature is greater than or equal to a temperature at a 
shallower depth, discards the corresponding depth. After 
proceeding through the entire interpolated trace in this manner, 
the trace is reduced to one in which there is only one depth 
corresponding to any given temperature, and it will be the 
shallowest depth at which this temperature appears. Although 
the program at this point has destroyed all information about 
inversions, this was not felt to be an unreasonable step. The 
raw data were all examined for inversions prior to the interpo- 
lation step, and in no case were inversions present greater than 
+0.1°C in five meters of depth. 

At this point, subroutine FINDZC is called, and it examines 
the temperature and depth information thus computed and chooses 
the depth of the first (shallowest) temperature which is less 
than or equal to ITC. It stores this depth as ZC, the depth of 
the center of the thermocline. This depth, along with all the 
remaining temperature and depth information, is returned to the 
main program. 

Using this information, the main program calls subroutine 
NONGA, which uses the procedure described in Section II to find 
the top and bottom of the thermocline through the use of 
subroutine GAUSS1 and GAUSS2, respectively. 
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For comparison with Tabata's results, the depth fluctuations 
of selected isotherms were also subjected to Fourier analysis. 
Depth information for this analysis was obtained in a manner 
similar to that described above, except that for depths of 
isotherms, subroutine FINDZS is called from INTERP after PURGE 
and the depths are determined in a manner similar to that used to 
find ZC . 

Fluctuations of isotherms and thermoclines with time are 
shown in Figures 3 and 4, respectively. 
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FLGURE 4 

Thermocline Depths 
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FIGURE 3 

Isotherm Depths 



V. COMPUTATION OF POWER SPECTRA 



The reworked data obtained as described in Section IV were 
subjected to Fourier analysis to produce a power spectrum. 
Subroutine RHARM, which is contained in the IBM Scientific 
Subroutine Library, was used for this analysis. This subroutine 
uses an algorithm, developed by Tukey and Cooley |l965 , to 
transform the depth fluctuations to a Fourier series of the form 

= N a n cos ncot+a n s.n ncot 

nro' 



where flCJ is the frequency of a given sinusoidal component in 
radians per hour, N is half the total number of data points plus 
one, t is time in hours, andT^is the instantaneous wave elevation 
above its mean position. The sum of the squares of the 
coefficients, a n +b n , is proportional to the power contained in 
the nth component, whose period is N/ (n) hours. ag is the mean 
value of the data, and to eliminate leakage of this large value 
into neighboring frequency components, all data was adjusted to 
a mean value of zero in subroutine AVERGE before RHARM was called, 
bg is always equal to zero, and is therefore ignored in subsequent 
manipulations . 

Once the Fourier coefficients were obtained, they were hanned 
by subroutine HANN. This subroutine used the following recursion 
formula to refine the coefficients: 

t 

a^ = 0.5a^ - 0.25a2 

a^ =-0.25a i _ 1 + O.Sa.^ - 0.25a^ + jL i = 2,3, ,N-1 

a N = “^*2^ a N-l + 0.5ajq 

The sine coefficients were similarly treated. Hanning the data is 



24 



a method of greatly reducing leakage of power from frequencies 
containing large amounts of power into neighboring frequencies. 
According to Bingham, et. al. fl967j if the Fourier transform is 
perfect for^J-^, CO 3,... and CO is none of the the 



-1 



leakage without harming is proportional to | CJj— M 

the coefficients reduces the leakage to a term proportional 

-3 



and harming 



to 



Cu r CO 



The hanned data was then plotted on a graph of power (spectral 
density) versus frequency. Although this presentation gave 
indications of the frequencies of the major components of the 
internal waves present, it was objectionably noisy. To reduce 
this noise, a method of power spectrum smoothing was utilized to 
produce the final spectrum. This method, originated by 
Jones Il965 was used in subroutine SMOOTH. The method consists 



l 

of first arranging the squares of the non-zero a and b coef- 
ficients in order and renaming them s, as follows: 

2 2 2 2 2 
a 0 * a l ,a l ,a 2 »"2 » 

S 1 * s 2 ,s 3 » s q ,s 5 > 

The smoothed spectrum is obtained by averaging these numbers in 
groups of 2m, except for the first group, which contains m, and 
the last group, which contains those numbers which do not fit 
into groups of (2p + l)m, where p is the number of groups of 
2m included in the data. 
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S - N/m 
0 

< 


s-j + 2 / s . 

^4 i=2 


S ± = N/m 


> 


rr> 1 


• 

• 




Sp = N/mp 


) S i 


(2p-T'm.1 



where rrip is the number of pairs remaining in the last step and 
N is the total number of Fourier coefficients used in the 
smoothing. By this method of smoothing, S Q will be centered at 
,25/N 8t, will be centered at 1/N bt, S 2 at 2/N 8t, and so 
forth. Care must be taken to choose m sufficiently large that 
significant smoothing takes place without choosing it so large 
that it eliminates detail from the final spectrum. 

Subroutine SMOOTH was run with m’s of 2, 3, 4, 5, 10, 20, 
and 50. Of these, the best results were obtained with m = 4, 
and this grouping factor was used in the final analysis. After 
the data were smoothed, subroutine PERIOD was called to generate 
counters to be used for the final graphic output. 

The data were broken into two groups, the first containing 
the first 512 data points (consecutive numbers 1 through 512) and 
the second containing the last 512 data points (consecutive 
numbers 139 through 651) . More confidence may be placed in 
spectral estimates if the spectral estimates are the averages 
of several sets of independent data. The nearest approximation 
to this which was possible with this data was the treatment which 
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was used. Had a sample size of 256 points been chosen, totally 
independent data sets could have been used, and confidence in 
the results would have been increased. The reduction in the 
number of data points by a factor of two, however, would have so 
significantly reduced the accuracy of the harmonic analysis that 
the increased confidence in the results would have been 
meaningless. For this reason, the use of two partially independ- 
ent data sets of 512 points each was decided upon. 

The two groups of hanned, smoothed data were then added 
together, term for term, by subroutine GROUP, and averaged. 

Plots of the results of this averaging for the top, center, and 
bottom of the thermocline are shown in Figures 5, 6, and 7. 

Analysis of the depth fluctuations of selected isotherms was 
carried on in the same manner as analysis of the thermocline 
data. The smoothed spectra of fluctuations of the isotherms 
analyzed are shown in Figures 8, 9, 10, 11, 12, and 13. For 
comparison, Tabata’s results for the 6°C and 10° C isotherms are 
shown in Figure 19. 

The complete computer program used in data reduction is 
given in the Computer Program section of this thesis. 
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FIGURE II 
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6°C&10°C ISOTHERM SPECTRA (tabata,1965). 
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VI. RESULTS 



In his paper, Tabata reports a strong semidiurnal component 
in the internal waves analyzed through depth fluctuations of the 
6 c C through 10° C isotherms, with greater power in higher 
frequencies noted in analysis of the deeper isotherms. This 
general behavior was also noted in the digitally analyzed data 
used in this thesis. 

A power peak centered at an 11.64 hour period was present 
in the spectra for all the isotherms analyzed, as well as in the 
spectra of the top, center, and bottom of the thermocline--the 
shallowest feature studied — also shows a peak corresponding to a 
period of 25.6 hours, probably associated with a principal lunar 
diurnal (0j_) tidal component. 

The power spectrum due to fluctuations of the bottom of the 
thermocline indicates a shift of power to higher frequencies, 
with peaks associated with periods of approximately 4.5 and 3.5 
hours, and a maximum energy associated with waves of 9.14 hour 
period. Tabata reports waves of periods equal to five hours 
and slightly less than four hours present for the deepest of 
the isotherms which he studied, with which these data agree . In 
no case, however, does he report waves of periods near nine 
hours. Relative peaks in spectral density having periods between 
nine and ten hours (most often at a period of 9.85 hours) were 
found associated with all features. Waves with periods of 
16.00 hours, which is very close to the inertial . period lof IS, 7 
hours at Ocean Station "P” were also noted in all spectra 
produced. 
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In the original plan for the analysis, the fluctuations of the 
thickness of the thermocline were also to be investigated, to 
determine whether the thermocline was oscillating in the funda- 
mental mode — in which the entire thermocline oscillates as a 
unit- -or mode two- -in which the top and bottom of the thermocline 
are 180° out of phase, resulting in a periodic thickening and 
thinning of the thermocline* s vertical extent. The spectra pro- 
duced by analysis of the thickness showed nearly constant spec- 
tral density throughout the entire range of frequencies studied, 
except for large contributions in the low frequency end of the 
spectrum. This suggested that the thermocline was oscillating in 
the fundamental mode, and the spectral densities at high 
frequencies were due to noise. To verify this, the correlations 
of the top, center, and bottom of the thermocline were computed. 

It was found that the depth fluctuations of each of these 
features with the others was greater than +.90 which confirms 
the fact that the entire thermocline was oscillating as a unit. 
This conclusion agrees with the results published by Tabata, who 
found phase differences on the order of - 20° between oscillations 
of the isotherms he analyzed. 

The spectrum associated with the top of the thermocline 
showed a very large contribution due to waves of very long 
period. Examination of Figure 3 shows the reason for this. The 
top of the thermocline was slowly sinking with time, presumably 
due to mixing in the upper layer. This effect, although small 
in magnitude , would have been identified in the analysis as an 
internal wave with a period equal to the total length of the 
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record, which in this case was 512 hours. In addition, during 
the period 16 August to 21 August, the water column, which was 
almost perfectly isothermal above and below the thermocline at 
al 1 other times, developed near-surface transient thermoclines . 

As the non-Gaussian analysis used to obtain thermocline depths 
did not eliminate transient thermoclines before computing the 
depths associated with the main thermocline, these near-surface 
transients resulted in the calculation of an abnormally large 
variance by subroutine GAUSS1, which had the effect of moving 
the calculated top of the thermocline closer to the surface than 
it actually was. The fact that the time period during which this 
occurred was included in both sets of 512 hourly observations 
meant that an additional large power component associated with 
a period of 512 hours was found in the analysis. This is an 
apparent effect, rather than a real one, and even when combined 
with the long-term sinking of the thermocline, fails to obscure 
the results in the frequency range of interest. 

Examination of Figure 3 also reveals that the center of the 
thermocline also showed a slight downward displacement with time, 
but of less magnitude than that of the top of the thermocline. 
This resulted in a relatively large amount of power associated 
with long-period internal waves at the center of the thermo- 
cline, but much smaller than the power contribution to the 
spectrum of the top of the thermocline . The bottom of the 
thermocline showed no change in mean depth with time, which 
indicates that surface layer mixing extended only to a depth 
slightly below the center of the thermocline. Although in most 
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cases large amounts of power were associated with low frequencies, 
these frequencies were not considered of interest, and are not 
indicated on the graphs of the final spectra. A complete presen- 
tation of all data produced by the analysis is presented in the 
COMPUTER OUTPUT section. 

Tabata's results for the 6°C and 10P C isotherms are presented 
in Figure 14 in a semilog plot for comparison with the results of 
the digital analysis. Although these curves differ in detail 
from those presented for the same isotherms analyzed digitally in 
this thesis, the overall form is strikingly similar. For both 
analyses of the 10°C isotherm, the greatest peak spectral density 
is associated with waves whose period is greater than the semi- 
diurnal tidal period. For both analyses of the 6°C isotherm, 
the spectral density is rather well- distributed throughout the 
higher frequency components, which also have a greater mean spec- 
tral density than the higher frequency components of the analyzed 
10' C isotherm data. This indicates that the higher frequency 
components showed an increase in amplitude with increasing depth, 
a behavior which was also noted by Haurwitz, Stommel, and Munk 
[l959] . 

The Brunt-Vaisala frequency associated with the thermal 
structure of this time series is 3.79 x 1Q - ^ sec - ’* - . This indi- 
cates that the thermal structure is capable of supporting inter- 
nal waves of any period greater than 26.4 seconds. The Nyquist 
frequency for a digital time series sampled at hourly intervals 
is 0.5 hr -1 . Thus, there is a large range of frequencies which 
may be present but of which the data does not permit investigation. 
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The causes of the internal waves with periods near 24, 16, 
or 12 hours seem obvious. The 24 and 12 hour periods are 
probably associated with internal waves which are being driven 
by diurnal or semidiurnal tidal forces. As mentioned previously, 
the 16 hour period seems to be associated with inertial lintarnall 
waves. The causes for the internal waves with periods between 
four and six hours — indicated both in this thesis and in Tabata’s 
study — are not so easily theorized. 

Another possible source of internal waves is baroclinic 
Rossby waves, but this type of wave is associated with extremely 
low frequencies, with periods on the order of days, so it would 
seem that these waves will not explain the higher frequency 
waves indicated by the analysis. Two other possible causes for 
these waves do exist, however. They might be long-period waves 
which have been reflected from a coast and have interacted with 
varying density gradients, which act as a weak filter, but this 
theory cannot be tested unless one knows the direction of 
propagation of the higher frequency internal waves. 

The other possible source of these waves is a non-linear 
interaction between surface pressure .disturbances , but in order 
to examine this as a possible cause, one must have records of 
both surface waves and weather. Although weather records exist 
for the time during which the bathythermographs were taken, no 
surface wave information is available . To be of value in 
investigating this process as a causative agent, the surface 
wave information would have to be very detailed, and wave records 
of the required detail can only be obtained in the deep ocean 
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through highly sophisticated and careful measurements which are 
not compatible with the mission of a weather ship. 

In conclusion, it seems that little can be said about these 
higher-frequency internal waves except that they exist. In order 
to say more, one would require a great deal more data than is 
available. Little has been done in the investigation of internal 
waves to date, and one can see that it represents a fertile 
subject for further research. 
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VII. RECOMMENDATIONS 



During the execution of this thesis, a number of changes in 
the format and amount of data were wished. Had these changes 
been made, the analysis of the data would have been simplified 
and the accuracy of the results enhanced. For these reasons, the 
following recommendations are made . 

In digitizing mechanical bathythermograph traces, it would 
greatly assist the user of the digitized data if the data 
through a thermocline were recorded as depths at equally spaced 
temperatures, instead of temperatures at equally spaced depths 
as is now done. This would provide an increased number of points 
in the region of maximum gradient, enabling the digitized data 
to be used to more accurately reconstruct the original trace. 

The method presently used to digitize expendable bathythermograph 
traces, which records a data point at the location where the 
trace departs from linearity by a specified amount, allows use 
of digitized data to almost perfectly reconstruct the trace. 

If salinity data were also available for analysis, it would 
be a simple task to compute depths of isopycnals using a 
computer. If terms to second order involving salinity and 
pressure were introduced to the already-known temperature 
information, depth-density profiles of much greater accuracy 
would result, and this in turn could be analyzed to produce 
internal wave information of higher quality. 

Closer spacing of observations would enable analysis of 
higher frequency components, since through aliasing the smallest 
wave period that can be detected using a sampling interval of one 
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hour is two hours. Also, as was suggested by Tabata, utilizing 
ilnta presently available over extremely long periods, with t as 
large as a year, would allow analysis of extremely low frequency 
internal waves. 

As was indicated in the Results section of this thesis, 
there is a need for information regarding the direction of 
propagation of the internal waves indicated. This could be got- 
ten using a cross array, consisting of two orthogonal arrays 

r i 

of the type described by Gilchrist 19661 for use in investi- 
gating the direction of propagation of surface waves. SeveraJ 
such arrays would be required to allow investigation of a large 
range of possible wavelengths, but it would seem that the 
information obtained as a result of such an expenditure of effort 
and equipment would be a highly significant contribution to the 
present knowledge of internal waves. 
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APPENDIX A 



EXAMPLE OF COMPUTATION 



The following four figures indicate graphically the operation 
of the subroutines used in computing the depths used in harmonica 
analysis. The data is from bathythermograph number 392 (Fig. 2) .• 

Figure 15 indicates the operation of subroutine INTERP in fit- 
ting straight line segments between raw data points, and the 
choosing of a temperature corresponding to the center of ther- 
mocline (ITC) . 

Figure 16 shows the operation of subroutine PURGE in elimi- 
nating small inversions and duplications of temperatures, reducing 
the data to a set of temperature -depth pairs from the surface to 
the bottom of the trace at an interval of 0.1°C. This figure 
also indicates the operation of subroutines FINDZC and FINDZS in 
determining the depth of the center of the non-Gaussian thermo- 
cline and the depths of the selected isotherms. 

Figures 17 and 18 show the thermal structures used in GAUSS1 
and GAUSS2 to obtain the depth of the top and bottom of the 
thermocline, respectively. 



46 



5 6 



TEMPERATURE (°C) 

7 8 9 10 11 12 



13 




FIGURE 15 
Subroutine I N T E R P 



47 



DEPTH (m) 



TEMPERATURE (°C) 

8 9 10 11 12 

T 




FIGURE 16 
Subroutine PURGE, 

Subroutine FINDZC, 
and 

Subroutine FINDZS 
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FIGURE 17 
Subroutine GAUSS1 
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FIGURE 18 
Subroutine GAUSS2 
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4.92 


35.10 


41.82 


30.37 


4.74 


30.51 


18.39 


24.97 


4.57 


13.55 


29.50 


46. 12 


4.41 


53.24 


40.96 


49.66 


4.27 


12.45 


17.95 


31.90 


4.13 


25.96 


42.56 


63.11 


4. CO 


23.16 


28.28 


41.16 


3.86 


10.73 


14.97 


28.30 


3.76 


10.97 


32.88 


45.66 


3.66 


17.16 


14.16 


41.44 


3.56 


20. 07 


31.26 


94.63 


3.46 


42 .61 


86.33 


101.53 


3.37 


20.63 


10.14 


13.16 


3.2e 


6 .05 


22.31 


6.96 


3.20 


12.41 


15.48 


34.15 


3.12 


4 . 19 


28.65 


26.31 


3.05 


9.05 


22.92 


27. C7 


2.98 


19.21 


10.70 


11.19 


2.91 


22.03 


26.70 


20. 16 


2.84 


11.30 


23.64 


68.37 


2. 7e 


13.07 


18.31 


12 . 14 


2.72 


20.97 


26.19 


24.85 


2.67 


8.59 


18.77 


18.18 


2.61 


15.35 


15.45 


8.85 


2.56 


18.28 


18.59 


8. e2 


2.51 


11.36 


2C.89 


23.45 


2.46 


3.21 


19.45 


24.60 


2.42 


18.38 


31.24 


47. 12 


2.37 


3 • 9C 


11.29 


4C. 13 


2.33 


5.74 


15.77 


8.63 


2.29 


12.89 


13.49 


15.45 


2.25 


20.05 


2 7.29 


27 .C4 


2.21 


13.89 


13.72 


23.94 


2.17 


14.80 


33.53 


42.12 


2.13 


6.94 


18.84 


49.99 


2.10 


17.86 


14.31 


6.35 


2.C6 


14.26 


37.02 


13.11 


2 . C 3 


6.93 


11.23 


3.57 


1.99 


5.42 


28.84 


15.93 



S3 



SPECTRAL ESTIMATES FCR GROUP NUMBER 2 



PER I CD 


TCP 


CENTER 


BOTTOM 


256. CC 


1873.16 


52.36 


29.24 


128.00 


1230.09 


34.40 


65.41 


64.00 


62 .73 


71.81 


39.17 


42.67 


25.24 


33.78 


19.51 


32.00 


7.52 


9.03 


20.11 


25. 6C 


90.39 


44.97 


9.19 


21.33 


59.56 


27.46 


39.47 


18.29 


43 .80 


33.56 


15.10 


16. CC 


15.79 


41.30 


41.12 


14.22 


4.C8 


9.56 


10.88 


12. 8C 


24.20 


14.31 


22.51 


11.64 


44.52 


58.22 


39.16 


10.67 


20.60 


28.23 


14.35 


9.85 


45.76 


83.83 


31.85 


9.14 


29.83 


78.72 


94.26 


8.53 


12 .49 


44.50 


30.68 


8. CO 


18.52 


36.38 


42.39 


7.53 


24. C4 


43.13 


22.91 


7.11 


31.31 


26.27 


48.27 


6.74 


10.81 


21.69 


18.50 


6.4C 


34.62 


26.15 


4.79 


6.10 


8.79 


5.32 


4.89 


5.82 


17.31 


61.58 


19.61 


5.57 


5.10 


13.72 


8.98 


5.33 


20.25 


27.41 


24.39 


5.12 


7.43 


28.91 


16.09 


4.92 


18.44 


24.20 


13.47 


4.74 


17. C9 


43.48 


26.44 


4.57 


22. C7 


22.05 


34.15 


4.41 


12.50 


25.19 


12.74 


4.27 


4.51 


16.64 


12.88 


4.13 


4.02 


6.62 


12. 08 


4.00 


8 .21 


20.35 


35.31 


3. ee 


9.86 


13.39 


34.24 


3.76 


8.19 


18. 1C 


4.96 


3.66 


3.97 


5.08 


30.45 


3.56 


3.49 


8.71 


23.32 


3.46 


9.54 


39.34 


40.22 


3.37 


2.46 


8.65 


8.75 


3.28 


8.32 


6.48 


8.51 


3.20 


7.52 


13.95 


11.41 


3.12 


10.27 


7.18 


16.97 


3.05 


18.60 


16.51 


13. 08 


2.98 


9. 50 


9.84 


18.43 


2.91 


5.25 


5.91 


5.05 


2.84 


5.94 


10.51 


9.47 


2.78 


8.10 


15.46 


18.04 


2.72 


6.22 


9.52 


15.05 


2.67 


12.62 


9.96 


9.56 


2.61 


2 .99 


9.36 


4.56 


2.56 


7.37 


6.54 


8.84 


2.51 


13.50 


12.56 


36.85 


2.46 


17.12 


13.96 


12.55 


2.42 


6.24 


12.82 


22.43 


2.37 


4.75 


13.57 


7.89 


2.33 


1.19 


20.42 


35.32 


2.29 


7. 40 


14. 15 


7.85 


2.25 


4.74 


9.45 


4.69 


2.21 


8.02 


8.11 


12.78 


2.17 


9.94 


7.15 


1.35 


2.13 


15.54 


7.77 


11.75 


2. 10 


18.48 


14.44 


11.56 


2.06 


0.76 


4.66 


13.29 


2.03 


0 .82 


4.41 


5.74 


1.99 


3.93 


3.27 


28.37 



54 



SPECTRAL ESTIMATES FOR GROUP NUMBER 2 



PER IOC 


11 DEG 


10 DEG 


9 DEC- 


256. CC 


52.16 


69.90 


108. 77 


128. OC 


26 .C 1 


31.23 


38.58 


64. CC 


54 .4C 


66.58 


75.30 


42.67 


27.93 


32.99 


36.67 


32. CC 


12.73 


13.47 


14.18 


25. 6C 


47.71 


44.67 


42.49 


21.33 


21. 5C 


21.45 


25.32 


18.29 


35.28 


32.62 


24.70 


16. CC 


40 .05 


45.20 


54 .94 


14.22 


8.69 


13.78 


18.46 


12. 8C 


13 .8e 


12 .51 


8.94 


11.64 


60.50 


53.13 


52.29 


10.67 


31.83 


29.07 


3C.83 


9.85 


74.73 


75.01 


68. C6 


9.14 


62.45 


62.82 


55.57 


8.53 


46.58 


43.35 


40.45 


8.00 


36.11 


33.86 


33.10 


7.53 


35.29 


35.05 


32.63 


7.11 


31 .CO 


26.05 


24.60 


6.74 


18.77 


20.93 


22.63 


6.4C 


29.29 


28.22 


22.97 


6. 1C 


3 .03 


4 . 18 


3.45 


5.82 


57.95 


55.82 


51.14 


5.57 


10.39 


12.71 


14.96 


5.33 


13.85 


16.17 


18.22 


5.12 


28.18 


27.08 


34.45 


4.92 


28.18 


28.90 


28.48 


4.74 


45.51 


49.04 


41.54 


4.57 


9.48 


11.09 


13.60 


4.41 


24 .oe 


25.39 


24.76 


4.27 


12.55 


14.16 


15.29 


4. 13 


10.37 


11.26 


1C. Cl 


4.00 


22.87 


21.64 


17.14 


3.88 


15.23 


13.69 


11.97 


3.76 


18.28 


19.52 


20.96 


3.66 


8.13 


7 .64 


5.53 


3.56 


9.43 


10.50 


13.65 


3.46 


19.05 


24.26 


32.10 


3.37 


9.19 


8.84 


9. 03 


3.28 


8.85 


7.8C 


8.40 


3.20 


14.66 


15.48 


14. 16 


3.12 


8.10 


7.28 


5.28 


3 . C 5 


14.92 


17.01 


16.74 


2.98 


9.45 


7.67 


8.27 


2.91 


11.37 


9.91 


11.06 


2.84 


6.44 


7.14 


10.95 


2.78 


13 .31 


13 .44 


12.51 


2.72 


9.73 


9.70 


9.37 


2.67 


9.18 


8.74 


9.78 


2.61 


12.80 


14.21 


16.39 


2.56 


2.42 


4.57 


8.24 


2.51 


8 .99 


9.65 


14.19 


2.46 


8 .45 


9.08 


11.27 


2.42 


11.11 


12.28 


13.77 


2.37 


5.68 


8.81 


9.51 


2.33 


13 .06 


15.17 


16 .70 


2.29 


9.79 


10.11 


11.41 


2.25 


3.36 


3.40 


7.66 


2.21 


13.41 


13.17 


13.00 


2.17 


2.21 


2.81 


6.88 


2.13 


4 .79 


4 .42 


5.49 


2. 1C 


12.13 


12.21 


9.89 


2.C6 


7.11 


8.91 


9.91 


2.03 


4 .67 


4.57 


4.76 


1.99 


4.76 


4.45 


3.54 



55 



SPECTRAL ESTIMATES FOR GROUP NUMBER 2 



PERIOD 


8 DEG 


7 DEG 


6 DEG 


256.00 


194.69 


213.73 


173.83 


128.00 


57.58 


129.74 


147.77 


64.00 


76 .76 


133.56 


193. 10 


42.67 


55.49 


100.73 


79.06 


32. OC 


21.23 


76.60 


65.44 


25.60 


42.05 


55.35 


18.11 


21.33 


21.95 


36.54 


40.63 


18.29 


11.65 


30.71 


25.92 


16. CO 


46.29 


66.71 


135.36 


14.22 


18.27 


19.70 


42.29 


12.80 


4.94 


6.52 


50.43 


11.64 


48. CC 


53.58 


48.18 


10.67 


38.98 


49.28 


38.86 


9.85 


63.58 


68.84 


55.95 


9.14 


55.29 


80.36 


57.35 


8.53 


39.34 


31.65 


7.71 


8.00 


33.78 


43.01 


60.06 


7.53 


30.41 


39.84 


28.44 


7.11 


17.66 


54.21 


116.32 


6.74 


23.74 


23.54 


15.47 


6.40 


20. 1C 


21.76 


36.51 


6.10 


9.99 


6.70 


4.73 


5. 82 


71.99 


61.25 


52.81 


5.57 


15.34 


36.78 


23.22 


5.33 


9.67 


6.59 


17.98 


5.12 


34.92 


35.62 


21.27 


4.92 


30.63 


28. 13 


28. 39 


4.74 


39.80 


23.08 


32.17 


4.57 


12.48 


47.87 


64.83 


4.41 


21.42 


65.73 


38.23 


4.27 


15.05 


26.20 


10.77 


4. 13 


6.34 


18.78 


14.18 


4.00 


15.61 


14.13 


20.43 


3. 8 8 


12.49 


15.41 


23.33 


3 • 7fc 


14.23 


30.86 


51.93 


3.66 


6.18 


16.89 


18.91 


3.56 


13.21 


33.64 


70.53 


3.46 


33.28 


42.80 


33.14 


3.37 


11.39 


11.50 


49.13 


3.28 


14.58 


15.97 


19.08 


3.20 


9.39 


21.56 


20.77 


3.12 


9.07 


25.88 


15.84 


3.05 


15.49 


22. 13 


19.59 


2.98 


10.36 


6.02 


9.61 


2.91 


11.80 


28.53 


43.65 


2.84 


7.61 


15.23 


18.90 


2.78 


8.62 


10.36 


20.08 


2.72 


11.20 


29.23 


23.79 


2.67 


18.19 


16.42 


13.42 


2.61 


11.65 


32.13 


48. Cl 


2.56 


7.39 


10.89 


22.35 


2.51 


17.48 


22.07 


56.28 


2.46 


10 .91 


28.46 


45.81 


2.42 


21.95 


48.08 


41.00 


2.37 


8.38 


7.11 


21.15 


2.33 


18.11 


35.27 


19.77 


2.29 


12.87 


10.35 


17.23 


2.25 


17.15 


9.62 


11.68 


2.21 


19.65 


12.33 


22.57 


2.17 


7.90 


18.77 


15.58 


2.13 


7.50 


10 .18 


35.26 


2. 10 


8.71 


4.04 


12.44 


2.06 


9.69 


14.64 


8.94 


2.03 


5.92 


9.51 


3.79 


1.99 


3.34 


19.49 


14.29 



56 



AVERAGEC SPECTRAL ESTIMATES 



NR 


PER ICC 


TCP 


CENTER 


8CTTCM 


0.5 


256.00 


1611.60 


113.37 


24.20 


1.0 


128. CO 


1166.87 


55.85 


54.93 


2 .0 


64.00 


91.42 


83.71 


27.22 


3.0 


42.67 


64.69 


29.58 


19.80 


4.0 


32. CC 


47.53 


20.40 


14.77 


5 . C 


25.60 


1 10.29 


50.83 


23.66 


6.C 


21.33 


80.94 


28.71 


48.73 


7 . C 


18.29 


44.66 


33.89 


32.28 


fl.C 


16. CC 


40.23 


64.42 


70.67 


9.0 


14.22 


26.10 


30.40 


22.52 


10.0 


12.80 


30.06 


26.59 


33.78 


11.0 


11.64 


96.46 


109 .93 


68.12 


12.0 


10.67 


33.56 


4 6.61 


28.24 


13.0 


9.85 


43.40 


93.51 


36.28 


14 .0 


9. 14 


27.96 


55 .63 


87.34 


15.0 


8.53 


16.37 


57.95 


26.39 


16.0 


8.00 


13 .02 


29 .93 


39.37 


17. C 


7.53 


29.98 


43 .93 


30.50 


18. C 


7.11 


32.42 


32.37 


45. 84 


1*5 ,0 


6 . 74 


12.29 


33.62 


11.73 


20. 0 


6 . 4 C 


47.12 


25.09 


5.10 


21.0 


6.10 


11.17 


6.83 


4.01 


22.0 


5.82 


11.60 


55.35 


23.84 


23. C 


5.57 


7.52 


18.66 


12.02 


24. C 


5.33 


19.73 


31.02 


17.48 


25. C 


5. 12 


11.21 


23.54 


10.10 


26.0 


4.92 


14.69 


22.72 


14.94 


27 .0 


4.74 


25.99 


44.32 


24.21 


28.0 


4.57 


22.62 


24.08 


43.84 


29.0 


4.4 1 


22.06 


42.40 


26.43 


30 .0 


4.27 


6.30 


15.57 


10.43 


31. C 


4. 13 


7.85 


14.23 


20.37 


32.0 


4 .00 


12.71 


23.95 


23.87 


33.0 


3.88 


11.18 


15.63 


29.13 


34.0 


3. 76 


10.30 


14.22 


17.35 


35.0 


3.66 


5.86 


6.54 


24.05 


36.0 


3. 56 


4.42 


12.81 


24.06 


37.0 


3.46 


10.33 


29.11 


43.42 


38.0 


3.37 


3.46 


13 .47 


8.18 


39.0 


3.28 


6.93 


8 .48 


8.11 


40 .0 


3.20 


11.78 


19.17 


16.63 


4 1.0 


3.12 


6.28 


7.76 


24.15 


42. C 


3. 05 


23.56 


11.92 


19.51 


4 3.0 


2.98 


15.56 


13.86 


1C. 94 


44.0 


2.91 


6.93 


13.24 


8.39 


45. C 


2. 84 


11.75 


8.46 


12.39 


46. C 


2.78 


7.45 


17.10 


22.30 


47. C 


2.72 


10.31 


12.39 


9.96 


48 .0 


2.67 


10.40 


12.31 


8.61 


49. C 


2.61 


5.14 


10 .94 


3.6 1 


50 .0 


2.56 


9.26 


12.73 


7.58 


51.0 


2.51 


11.91 


10.45 


22. 76 


52. C 


2.46 


15.58 


10 .50 


12. C8 


53 .0 


2.42 


6.92 


13.79 


18.09 


54.0 


2.37 


3.95 


7.77 


9.58 


55. C 


2.33 


1.79 


13.99 


21.94 


56. C 


2.29 


6.64 


15.37 


14.00 


57.0 


2.25 


9.22 


12 .06 


6.94 


58. C 


2.21 


7.71 


10.69 


10.59 


59. C 


2.17 


9.89 


9.37 


6. 80 


60. C 


2.13 


15.85 


8.97 


10.22 


61. C 


2. 1C 


14.94 


15.86 


15.74 


62 .0 


2.06 


1.21 


5.88 


12.63 


63. C 


2.03 


1.48 


5.70 


9.28 


64.3 


1.99 


3.29 


2.45 


26.72 



57 



AVERAGEC SPECTRAL ESTIMATES 



NR 


PERICO 


11 DEG 


1C DEG 


9 DEG 


0.5 


256. OC 


162.9 1 


165.76 


159.68 


1.0 


128. OC 


AA.92 


51.19 


65.29 


2.0 


6 A • 00 


67.77 


82.38 


94.57 


3.0 


42.67 


21.71 


25 .75 


29.99 


A .0 


32.00 


19.92 


21.37 


22.67 


5.0 


25.60 


56.58 


5A .65 


52.24 


fc.C 


21.33 


25.39 


23.70 


23.31 


7.C 


18.29 


AO. 27 


37.08 


31.31 


8.0 


16.00 


69.02 


69.30 


76.93 


9.0 


1A.22 


31.83 


34.73 


38.67 


10.0 


12.80 


25.37 


26.21 


25.85 


11. c 


11. 6A 


115.22 


107.12 


104.43 


12.0 


10.67 


A5.8C 


43.77 


45.55 


13.0 


9.85 


80.59 


82.38 


8C.97 


14. C 


9. 14 


A7.02 


A7.02 


42.24 


15.0 


8.53 


56.39 


5A.23 


53.42 


16.0 


8.00 


31.55 


28. AO 


27.25 


17.0 


7.53 


39. 2A 


38.39 


36.31 


18.0 


7.11 


33.06 


28.83 


26.87 


19.0 


6.7A 


35.61 


37.90 


40.03 


20.0 


6. AC 


2A. 61 


24.21 


20.84 


21.0 


6.10 


6.70 


7 . 1 A 


5.89 


22.0 


5.82 


53.21 


51 .83 


46.86 


23.0 


5.57 


15.05 


16.19 


16.22 


24.0 


5.33 


17.53 


20.50 


23.70 


25. C 


5.12 


20.67 


20. A7 


25.93 


26.0 


A. 92 


26.93 


27.85 


28.97 


27.0 


A. 74 


42.42 


43.27 


37.00 


28. C 


A. 57 


12.75 


14.84 


16.32 


29.0 


A • A 1 


35.58 


38.27 


38.00 


30. C 


A. 27 


12.98 


13.80 


14.15 


31.0 


A. 13 


17.97 


17.89 


16.36 


32.0 


A. 00 


22.9 1 


22.90 


20.63 


33. C 


3.88 


15.59 


14.99 


13.49 


34.0 


3.76 


13. 9A 


14.86 


15.79 


35.0 


3.66 


9.60 


9.15 


7.98 


36. C 


3.56 


13.33 


14 .41 


17.13 


37.0 


3.A6 


18. A6 


22.28 


25.80 


38.0 


3.37 


16.06 


16.02 


13.46 


39. C 


3.28 


10. A 3 


9.69 


8.32 


AO .0 


3.20 


15.16 


16.00 


14.42 


A 1 .0 


3.12 


6.59 


6.16 


4.83 


A2.0 


3.05 


10.39 


10.86 


11.14 


A3 .0 


2.98 


13.66 


13.17 


14.02 


A A . 0 


2.91 


15.85 


14.82 


14.36 


A5.0 


2. 8A 


8.07 


7.76 


9.75 


A6.0 


2.78 


14.55 


15.12 


13.94 


A 7 . C 


2.72 


11.78 


11.52 


12.65 


A8.0 


2.67 


10.37 


10.05 


10.09 


A9.0 


2.61 


13.77 


13.65 


14.83 


50.0 


2.56 


7.85 


10.26 


14.12 


51.0 


2.51 


8.62 


8.10 


10.24 


52.0 


2.A6 


10.28 


9.61 


10.20 


53.0 


2.A2 


13.99 


13.90 


14.68 


5A.0 


2.37 


A. 35 


5.72 


5.93 


55. C 


2.33 


10.87 


11.29 


11.48 


56. C 


2.29 


8.66 


9.62 


11.64 


57.0 


2.25 


6.7A 


7.57 


11.46 


58.0 


2.21 


13.78 


13.64 


13.43 


59.0 


2.17 


A. 88 


5.52 


9.06 


60.0 


2.13 


5.27 


5.15 


6.04 


61.0 


2. 1C 


11.61 


11.83 


11.01 


62.0 


2.06 


5.02 


6.65 


8.12 


63.0 


2.03 


A • 6 1 


A .95 


4.85 


64.3 


1.99 


3. AO 


3.35 


3.83 
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AV ER AG E C SPECTRAL ESTIMATES 



NR 


PER ICC 


8 DEG 


7 DEG 


6 DEG 


C . 5 


256. CC 


179.73 


196.59 


2C4 .59 


l.C 


128. 0C 


93.43 


177.30 


2C9.18 


2 .0 


64.00 


95.89 


148. CO 


182.67 


3.0 


42.67 


44.80 


74 .02 


91.57 


4.C 


32. OC 


32.45 


65.40 


64.49 


5.C 


25.60 


52.30 


68.61 


48.65 


6.C 


21.33 


14.06 


47.53 


44.27 


7 . C 


18.29 


20.34 


68.38 


52.25 


8 • C 


16. OC 


76 .66 


83.28 


124.55 


9.C 


14.22 


48.05 


45.58 


53.42 


10. C 


12.80 


28.13 


45.77 


91. 1C 


1 l.C 


11.64 


99.13 


1C8.10 


1C8.25 


12. C 


10.67 


54.81 


65.62 


49.64 


13.0 


9.85 


73. 9C 


80.53 


5e.43 


14. C 


5 . 14 


41.56 


69 .44 


53.38 


15.0 


8. 53 


58.27 


41.26 


14.81 


16.0 


8.00 


30.92 


42.27 


67.63 


17. C 


7.53 


32.71 


43.83 


54.71 


18. C 


7.11 


20.83 


46.77 


1C1.24 


19.0 


6.74 


40.58 


33 . 14 


22.06 


2 C.C 


6.4C 


20.88 


16.98 


28.16 


21.0 


6.10 


7.49 


7.14 


5.53 


22. C 


5.82 


66 .69 


82.23 


73.71 


23 .0 


5.57 


16.37 


36.31 


3C.40 


24.0 


5.32 


13.82 


22.67 


39.86 


25. C 


5.12 


28.61 


22 .61 


13.95 


26. C 


4.92 


32.36 


34.97 


29.38 


27.0 


4.74 


35.15 


20 .73 


28.57 


28. C 


4.57 


13.02 


38.68 


55.47 


29.0 


4.41 


37.33 


53.35 


43.94 


30.0 


4.27 


13.75 


22.07 


21.33 


31. C 


4.13 


16.15 


30.67 


38.65 


32.0 


4.CC 


19.39 


21.20 


3C.8C 


33. C 


3.88 


11.61 


15.19 


25.82 


34. C 


3.76 


12.60 


31.87 


48.79 


35.0 


3.66 


11.67 


15.52 


3C. 17 


36. C 


3.56 


16.64 


32 .46 


82.58 


27.0 


3.46 


37.94 


64.57 


67.34 


38.0 


3.37 


16 .01 


10.82 


31.14 


39. C 


3.28 


10.31 


19.14 


13.02 


4C.C 


3.20 


10.90 


18.52 


27.46 


41.0 


3.12 


6.63 


27.26 


21.08 


42. C 


3.05 


12.27 


22.53 


23.33 


43.0 


2.98 


14.78 


8.36 


1C.4C 


44. C 


2.91 


16.91 


27.62 


31.91 


45. C 


2.84 


9.46 


19.44 


43.64 


46. C 


2.78 


10.85 


14.34 


16. 1 1 


47. C 


2.72 


16.09 


27.71 


24.32 


48. C 


2.67 


13.39 


17.59 


15. 8C 


49. C 


2.61 


13.50 


23.79 


28.43 


50. C 


2.56 


12.83 


14.74 


15.59 


51.0 


2.51 


14.42 


21.48 


39.87 


52.0 


2.46 


7.06 


23 .96 


35.20 


53.0 


2.42 


20. 16 


39.66 


44.06 


54. C 


2.37 


6.14 


9.20 


30.64 


55. C 


2.33 


11.93 


25 .52 


14.20 


56. C 


2.29 


12.88 


11.92 


16.34 


57. C 


2.25 


18.60 


18.45 


19.36 


58. C 


2.21 


16.77 


13.02 


23.25 


59. C 


2.17 


11.35 


26. 15 


28.85 


60. C 


2.13 


7.22 


14.51 


42.62 


61 .C 


2. 1C 


13.29 


9.18 


9.39 


62.0 


2.06 


11.97 


25.83 


11.03 


63. C 


2.03 


6.42 


10 .37 


3.68 


64.3 


1.99 


4.38 


24.17 


15.11 
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MAIN PROGRAM 

THIS PROGRAM READS TEMPERATURES AT EVEN DEPTH INCREMENTS 
OF FIVE METERS, CALLS SUBROUTINES TO FIND DEPTHS OF THE 
TOP, CENTER, AND BOTTOM OF THE THERMOCLINE, AND THE 
DEPTH OF SIX SELECTED ISOTHERMS. THESE DEPTHS ARE SUB- 
JECTED TO A FAST FOURIER ANALYSIS TO DETERMINE THEIR 
FOURIER COEFFICIENTS, WHICH ARE THEN HANNED, SMOOTHED IN 
FREQUENCY BANDS, AND STORED FOR EACH GROUP OF 2**N DATA 
CARDS. THE PROGRAM MAY ALSO BE USED TO PRODUCE GRAPHIC 
DISPLAYS OF FLUCTUATIONS OF DEPTHS WITH TIME, AND THE 
FINAL AVERAGED POWER SPECTRA. 

VARIABLES USED IN THIS PROGRAM ARE — 

I GROUP— THE NUMBER OF GROUPS OF 2**N DATA CARDS USED 
IN EACH ANALYSIS 

IT— THE ARRAY OF TEMPERATURES IN EACH TRACE AS READ 
FROM DATA CARDS 

Z11,Z10,Z09,Z08,Z07,Z06 — THE DEPTHS OF THE SIX 
SELECTED ISOTHERMS, AS CHOSEN IN SUBROUTINE INTERP 

KN — THE SUBSCRIPT OF THE DEPTH CHOSEN AS THE INITIAL 
APPROXIMATION OF THE CENTER OF THE THERMOCLINE 

ZT — THE DEPTH OF THE TOP OF THE THERMOCLINE 

ZC — THE DEPTH OF THE CENTER OF THE THERMOCLINE 

ZBT — THE DEPTH OF THE BOTTOM OF THE THERMOCLINE 

J— USED AS THE SUBSCRIPT OF THE DEPTHS OBTAINED FOR 
EACH DATA CARD 

M — THE GROUPING NUMBER, A MEASURE OF THE NUMBER OF 
INDIVIDUAL FOURIER COEFFICIENTS INCLUDED IN EACH 
SPECTRAL BANDWIDTH. NOTE THAT M MUST BE SPECIFIED IN 
THE MAIN PROGRAM 

A J , B J , Z J — COUNTERS ASSOCIATED WITH THE CONSECUTIVE 
NUMBER OF THE SPECTRAL DENSITIES, FOR USE IN GRAPHICAL 
OUTPUT 

END1 1 , ENDIO , ENDT , ENDC , ETC . — FINAL SPECTRAL ESTIMATES 
USED TO OBTAIN AN AVERAGE POWER SPECTRUM 

SUBROUTINES CALLED FROM MAIN PROGRAM ARE 
INTERP 
NONGA 
AVERGE 
HANN 
SMOOTH 
PERIOD 
GROUP 

RHARM — FOR DETAILS CONCERNING SUBROUTINE RHARM SEE THE 
DOCUMENTATION FOR THE IBM SCIENTIFIC SUBROUTINE 
LIBRARY 

DIMENSION IT( 200) ,TX( 200) ,X(200) . APER ( 200 ) , A J ( 200 ) , 

*ZJ (200 ),BJ (200 >,Z 11(600 ) , Z 10( 600 > , Z09( 600 ) , Z08( 600) , 
♦Z07I600) , Z06 (600 ) ,B11 (600 ) , BIO ( 600 ) , B09 ( 600) , B08 ( 600 ) , 
*BC7( 600), 806(600) ,ZT( 600) ,ZC ( 600 ) ,ZBT ( 600 ) ,AZT(600), 
*AZC(600) , AZBT (600 ) , AZ IN (600 ) , AZOUT ( 600 ) , Z I N( 600) , 

♦BOUT (600) ,ARGIN(600) , ARGOUT (600 ) , INV ( 64) , S <64 ) . 

* ENDT ( 200 ) , ENDC ( 200 ), ENDS ( 200 ), ENDI 1(200) , ENDI 0(200), 

♦ END09 { 200 ) , END08 ( 200 ) , END07( 200 ) , END06( 200 I , T( 200 ) , 
*Z(200) 

REAL L ABL 1/ * 11 */,LABL2/* 10 */,LABL3/* 9 •/, 

♦LABL4/' 8 ' / , L ABL 5/ ' 7 */,LABL6/» 6 ' / , LA8L7/ ' TOP ' 
*/,LABL8/'MIDL'/,LABL9/*B0TM' / 

REAL*8 ITIK12)/' BOX 8009 AVERAGED POWER SPECTRUM OF 
♦FOURIER COEFFICIENTS DATA INPUT HAS MEAN OF ZERO H, 



*S4,ZT'/ 

RE AL*8 I T I 2 ( 12 )/* BOX 80D9 AVERAGED 
♦FOURIER COEFFICIENTS OATA INPUT 

*S4 ,ZC • / 

RE AL*8 I T I 3 I 12)/* BOX 8009 AVERAGED 
♦FOURIER COEFFICIENTS DATA INPUT 

*S4 , Z B * / 

REAL* 8 ITITL1II2)/* BOX 
♦TIONS OF THE 11 DEGREE 
♦NALYS IS'/ 

REAL* 8 I T I TL2 ( 12 I / • BOX 
♦TIONS OF THE 10 DEGREE 
♦NALYS IS'/ 

R EAL* 8 ITITL3I 12) /'BOX 
♦TIONS OF THE 9 DEGREE 
♦NALYS I S • / 

REAL* 8 ITITL4! 12) /'BOX 
♦TIONS OF THE 8 DEGREE 
♦NALYSI S'/ 

REAL*8 ITITL5I 12) /'BOX 
*T IONS OF THE 7 DEGREE 
♦NALYSI S'/ 

R EAL * 8 ITITL6I 12 ) /'BOX 
♦TIONS OF THE 6 DEGREE 
♦NALYSI S'/ 

DO 3 I GROUP= 1 » 2 

DO 1 J= 1 , 5 12 

READ ( 5 , 2 ) I I T ( I ) , 

2 FORMAT ( 13X, 2013) 

CALL INTEPPIIT,20 
♦Z06, J ) 

CALL NONGA ! Z, T » ICO UN T , 

1 CONTINUE 

CALL AVERGE 
AVERGE 
AVERGE 
AVERGE 
AVERGE 
AVERGE 
AVERGE 
AVERGE 
AVERGE 



POWER SPECTRUM OF 
HAS MEAN OF ZERO H, 

POWER SPECTRUM OF 



80 FOURIER ANALYSE 
CENTIGRADE ISOTHERM 

80 FOURIER ANALYSI! 
CENTIGRADE ISOTHERM 

80 FOURIER ANALYSI! 
CENTIGRADE ISOTHERM 

80 FOURIER ANALYSI! 
CENTIGRADE ISOTHERM 

80 FOURIER ANALYSI! 
CENTIGRADE ISOTHERM 

80 FOURIER ANALYSI! 
CENTIGRADE ISOTHERM 



1 = 1 , 20 ) 

T , Z , I COUNT , KN, Zii,Z10,Z09,ZC8,Z07, 
ZT,ZC,Z8T,KN,J) 



OF 


ZERO H, 


OF 


FLUCTUA 
RHARM A 


OF 


FLUCTUA 
RHARM A 


OF 


FLUCTUA 
RHARM A 


OF 


FLUCTUA 
RHARM A 


OF 


FLUCTUA 
RHARM A 


OF 


FLUCTUA 
RHARM A 



( Z T , AZ T ) 

CALL AVERGE (ZC.AZC) 

CALL AVERGE ( ZBT , AZBT ) 

CALL AVERGE (Zil,Bli) 

CALL AVERGE (Z10,B10) 

CALL AVERGE (Z09,B09) 

CALL AVERGE <Z08,B08) 

CALL AVERGE (Z07,B07) 

CALL AVERGE IZ06,B06) 

CALL RHARMIAZT,8, INV, S, IFERR) 
CALL RHARM(AZC,8, INV, S, IFERR) 
CALL RHARM(AZBT,8 , INV, S, IFERR) 
CALL RHARM(B11,8, INV, S, IFERR) 
CALL RHARM(B10,8, INV, S, IFERR) 
CALL RHARM(B09,8, INV, S, IFERR) 
CALL RHARMIB08 ,8, INV, S, IFERR) 
CALL RHARM(B07,8, INV,S, IFERR) 
CALL RHARM(B06,8, INV, S, IFERR) 
CALL HANN (AZT,ZT) 

CALL HANN (AZC,ZC) 

CALL HANNI AZBT, ZBT) 

CALL HANN (Bll,Zil) 

CALL HANN <B10,Z10) 

CALL HANN (B09,Z09) 

CALL HANN <B08,Z08) 

CALL HANN (B07,Z07) 

CALL HANN (B06,Z06) 

M=4 



CALL SMOOTH ( ZT , M , I TE S T , AZ T , A J , I PTS ) 
CALL SMOOTH ( ZC , M , I T EST , AZC , A J , I PTS ) 
CALL SMOOTH ( ZBT , M , I T EST , AZBT , A J , I PT S ) 
CALL SMOOTH! Z 11, M, I TEST, B 11, A J, I PTS) 
CALL SMOOTH! Z10,M, I TEST, BIO, AJ , I PTS > 
CALL SMOOTH! Z09 , M , I TE ST , B09 , A J , I PTS ) 
CALL SMOOTH! Z08,M,ITEST,B08,AJ,IPTS) 
CALL SMOOTH! ZO 7 , M , I T EST , B07 , A J , I PTS) 
CALL SMOOTH! Z06 , M , I TE ST , B06 , A J , I PT S ) 
CALL PERIOD!M, I T E ST , A PER , A J , B J ) 

WRITE !6, 30 ) 
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30 FORMAT ( 1H 1 ) 

WRITE(6»6) I GROUP 

6 FORMAT (//32X, • SPECTRAL ESTIMATES FOR GROUP NUMBER 
*11,//) 

WRITE (6,7) 

7 F0RMAT(30X,' PERIOD* ,10X,* TOP* , AX , ' CENTER* , AX , 

♦•BOTTOM* ,/ ) 

DO 8 1=1 , IPTS 

WRITE (6,9) APER(I),A2T(I),AZC(I ) ,AZBT(I ) 

9 F0RMAT(3CX»F6.2»AX,3F10.2) 

8 CONTINUE 

WR ITE ( 6, 30 ) 

WR IT E (6,6 ) I GROUP 
WRITE(6,10) 

10 F0RMATI30X, 'PERIOD*, 7X, • 11 DEG*, AX, *10 DEG*,5X,«9 DEG* 
* , / ) 

DO 11 1=1, IPTS 

WRITE! 6,9 ) APER( I ) , B 1 1( I ) ,B10( I ) , B09 ( I ) 

11 CONTINUE 

WR I TE ( 6, 30 ) 

WRITE ( 6, 6 ) I GROUP 
WRITE ( 6, 13 ) 

13 FORMATI30X, 'PERIOD* ,8X,* 8 DEG*,5X,*7 DEG*,5X,*6 DEG*, 
*/ ) 

DC 1A 1=1, IPTS 

1A WRITE(6,9) APER( I ) ,B08( I ) ,B07( I) ,B06( I ) 

DO 26 1=1, IPTS 

ENDT ( ( IGROUP-1 )*IPTS+ I )=AZT< I ) 

ENDC ( ( I GROUP- 1)*IPTS+I)=AZC( I ) 

ENDB ( ( IGROUP-1 )*IPTS+I )=AZBT( I ) 

END11( ( IGROUP-1)* IPTS+I)=B11( I) 

END10 ( ( IGROUP-1 )*IPTS+ I )=B10( I ) 

E ND09 ( ( IGRCUP-1)*IPTS+I) = B09( I) 

E NDO 8 ( ( IGROUP-l)*IPTS + I ) = B08( I ) 

END07 ( ( IGROUP-1 ) * IPTS+ I )=B07( I ) 

E ND06 ( ( IGRCUP-1)*IPTS + I) = B06( I ) 

26 CONTINUE 
3 CONTINUE 

CALL GROUP(AJ,BJ, APER, IPTS, ENDT, ENDC, ENDB, ENO 11, END10, 
*END09 » END08, END07 , END06 ) 

CALL DRAW FOR ISOTHERMS 
CALL DRAW THERMOCLINES 

INCREASE SIZE OF ELEMENTS OF ENDT, ENDC, AND ENDB 
FOR GRAPHING 
DO 999 1 = 1, IPTS 
ENDT ( I ) = 1 . 5*ENDT ( I ) 

ENDC (I )=1 . 5*ENDC ( I) 

999 EN0B( I )=ENDB( I )*1 .5 
WRITE(6»30) 

CALL DRAW ( IPTS ,B J ,END11 ,0 ,0,LABL1 , IT ITL1 , 1 . ,20. ,0 ,0, 2 , 
*2 , 5, 8 , 0, I L ) 

CALL DRAW(IPTS,BJ,END10,0»0»LABL2, ITITL2, 1., 20. ,0,0, 2, 
*2,5,8, 0»IL) 

CALL DRAW ( I PTS , B J , END09 , 0 , 0, L ABL 3 , I T I TL3 , 1 . , 20. ,0,0,2, 

*2 , 5 , 8 , 0, 1 L ) 

CALL DRAW ( I PTS , B J , END08 ,0 , 0 , LABLA , I T I TLA , 1 . , 20 . ,0 , 0 , 2 , 
*2 ,5,8,0, IL) 

CALL DRAW( I PT S , B J , END07 ,0 , 0 , L ABL5 , ITITL5, 1 . , 20 . , 0 , 0 , 2 , 
*2 ,5, 8, 0, I L ) 

CALL DRAW ( IPTS , BJ , END06 ,0 , 0, L ABL6 , IT I TL6 , 1 . , 20. ,0 , 0, 2 , 
*2,5, 8,0, IL) 

CALL DRAW! IPT S , BJ , ENDT, 0, 0, LABL7 , I T 1 1 , 1 . , 20. , 0, 0, 2, 2, 
*5 ,8,0, IL) 

CALL DRAW ( IPT S , B J , ENDC , 0 , 0 , LABL8 , 1 T 1 2 , 1 . , 2 0. ,0,0, 2, 2, 
*5 ,8,0, IL) 

CALL DRAW! I PTS , B J , ENDB,0 , 0 , LABL9 , IT 13 , 1 . , 20 . , 0 ,0 , 2, 2 , 
*5,8,0, IL) 

STOP 

END 
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SUBROUTINE INTERPI IT, I PTS , TX , X , ICOUNT , KN , Zll ,Z10, Z09, 
*Z08,ZC7,Z06,KK) 

THIS SUBROUTINE, CALLEO FROM THE MAIN PROGRAM, IS USED 
TO LINEARLY INTERPOLATE BETWEEN OATA POINTS GIVEN AT 
FIVE METER INCREMENTS TO DEPTHS CORRESPONDING TO TEMP. 
INCREMENTS OF 0.1 DEGREE CENTIGRADE. THIS PROGRAM DOES 
NOT CONSIDER SMALL NEAR-SURFACE TRANSIENTS OR INVERSIONS 

VARIABLES USED IN THIS SUBROUTINE 

IT — INPUT TEMPERATURES, IN TENTHS OF A DEGREE, FROM 
MAIN PROGRAM 

TX — THE ARRAY OF INTERPOLATED TEMPERATURES 

X — THE ARRAY OF INTERPOLATED DEPTHS ASSOCIATED WITH TX 

ICCUNT — THE TOTAL NUMBER OF ELEMENTS IN X AND IN TX 

I PTS — THE NUMBER OF TEMPERATURES PER IRACE AS READ 
FROM DATA CARDS. IN THE CASE OF THIS DATA IPTS=20 

IT1, IT2, IT3, IT4, IT5, I T6 — THE TEMPERATURES, IN TENTHS 
OF A DEGREE, OF THE ISOTHERMS WHOSE DEPTH FLUCTUATIONS 
ARE OF INTEREST. THESE VALUES MUST BE SPECIFIED IN 
THIS SUBROUTINE, AND ARE ASSOCIATED WITH Z11,Z10,Z09, 
Z08, Z07 , Z06 RESPECTIVELY IN LATER COMPUTATIONS. 



T , TT , I TT — ARRAY NAMES ASSIGNED TO 
TEMPERATURES AT VARIOUS POINTS IN 



THE INTERPOLATED 
THE COMPUTATIONS 



Z,ZZ — ARRAY NAMES ASSIGNED TO THE INTERPOLATED DEPTHS 

ITO — AN ARRAY USED TO FIND A TEMPORARY TOP, BOTTOM, 
AND CENTER OF THE THERMOCLINE. 

ITTOP — THE TEMPERATURE OF THE TEMPORARY TOP OF THE 
THERMOCLINE. DEFINED IN THIS SUBPROGRAM AS THE DEPTH 
BELOW WHICH THE TEMPERATURE FIRST DECREASES MORE THAN 
0.3 DEGREES IN FIVE METERS OF DEPTH. 

ITBOT — THE TEMPERATURE OF THE TEMPORARY BOTTOM OF THE 
THERMOCLINE. DEFINED AS THE FIRST TEMPERATURE LESS 
T FAN ITTOP ABOVE WHICH THE TEMPERATURE DECREASES LESS 
THAN 0.2 DEGREES IN FIVE METERS OF DEPTH 

ITC — THE ARITHMETIC MEAN OF ITTOP AND ITBOT, USED IN 
SUBROUTINE FINDZC TO DETERMINE THE SUBSCRIPT OF THE 
DEPTH VALUE CORRESPONDING TO THE ESTIMATED CENTER OF 
ThE THERMOCLINE 



SUBROUTINES 
PURGE 
FINDZC 
F INOZS 



CALLED FROM THIS SUBROUTINE ARE 



999 



DIME NS I ON IT(1),TX(1),X(1) , Zl 1 ( 1 ) , Z10( 1 ) ,Z09( 1 ) ,Z08( 1 ) , 
*Z07( 1 ) , Z06 ( 1 ) ,T( 200) , TT( 200) , Z ( 200 ) , Z Z ( 200 ) , I TT ( 2 00 ) , 

* I TO { 200 ) 

ILESS=IPTS-I 
GENERATE Z ( I ) 

DC 999 I = 1 , I PT S 
Z( I ) =5 . * ( I -1 ) 

CONTINUE 
I E PS=0 
IT1-U0 
I T 2= 100 
I T3=90 
I T4= 80 
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IT 5= 70 
I T 6=60 

DO 90 1 1 = 1 » I LESS 
I TO ( II )=IT ( I I+1)-IT( I I ) 

IFCITOIID.LT. -31 GO TO 91 

90 CONTINUE 

91 I T TO P= IT C 1 1 ) 

DO 92 JJ=10» I LESS 
I T0( J J )=IT(JJ+1)-IT(JJ) 

I F ( I TO ( JJ) .GT.-2) GO TO 93 

92 CONTINUE 

93 I TBOT= I T ( J J ) 

ITC=< ITTOP+ITBOTI/2 

EXAMINE I T ( I ) VERSUS ITU + l). IF SHALLOWER TEMP. 

IS GREATER THAN NEXT TEMPERATURE BELOW IT IN THE 
WATER COLUMN, INTERPOLATE. IF NOT, STORE. 

DO 1000 1=1, I PTS 

J=I+IEPS 

K=0 

I F ( I .EQ.IPTS) GO TO 71 
IFCITCD .GT.ITC DID GO TO 5 
I TT( J 1 = I T C I ) 

ZZ(J) = Z( I ) 

GO TO 1000 

5 I TDI F= IT ( I ) — IT ( 1 + 1 ) 

ZDIF=5 •/ ITDIF 
MAX= ITDIF+J-1 
DO 1001 J J=J , MAX 
ZZ(JJ)=( JJ-J )*ZDIF+Z( I) 

ITT( JJ)=IT(I)-K 
K=K+ 1 

1001 CONTINUE 

I EPS=MAX- I 
1000 CONTINUE 

INT= ILESS+J 
71 ITT(J)=IT( I) 

ZZ(J)=Z( I) 

100 DO 8 1=1, J 
TT ( I ) = ITT ( I ) 

8 CONTINUE 
I COUNT =J-1 

CALL PURGE ( TT , ZZ, ICOUNT , TX , X ) 

CALL FINDZC(TX,ITC,KN, ICOUNT) 

CALL FINDZSITX.X, IT1, IT2, IT3, IT4, IT5, IT6, Z 1 1 , Z 10 , Z09, 
*Z08,Z07,Z06,KK) 

20 RETURN 
END 



SUBROUTINE PURGE ( TT, ZZ , ICOUNT , TX , X) 

SUBROUTINE PURGE TAKES AS INPUT THE ARRAY OF TEMPERATURE 
AND DEPTH INFORMATION FROM SUBROUTINE INTERP AND 
EXAMINES IT FOR DUPLICATION OF TEMPERATURE, CONSIDERING 
A POSITIVE TEMPERATURE-DEPTH GRADIENT TO BE ISOTHERMAL. 
(THIS WAS JUSTIFIED BY THE PARTICULAR DATA INVOLVED IN 
THE WRITING OF THIS SUBPROGRAM.) THE SUBROUTINE STORES 
EACH TEMPERATURE IN THE TRACE AT THE SHALLOWEST DEPTH 
AT WHICH IT APPEARS. 

VARIABLES USED IN THIS SUBPROGRAM 

TT , TX , ZZ , X — AS DEFINED IN SUBROUTINE INTERP 

ICOUNT — AS OUTPUT, THE NUMBER OF TEMPERATURES OR DEPTH 
PAIRS REMAINING AFTER DUPLICATIONS HAVE BEEN REMOVED. 

DIMENSION TT ( 1 ) , TX ( 1 ) , ZZ ( 1 ) , X ( 1 ) 

TX( l)=TT( 1 ) 
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DO 1 1=1, ICOUNT 

IFCTTt I+ll .GE.TT( I ) ) GO TO 1 

TX( J)=TT( 1+1) 

X(J)=ZZ( 1+1) 

J = J+1 

1 CONTINUE 

ICOUNT = J-2 

RETURN 

END 



SUBROUTINE FINDZS(TX,X,IT1,IT2,IT3,IT4,IT5,IT6,Z11, 
*Z10,Z09,Z08,Z07,Z06,KK) 

THIS SUBROUTINE TAKES THE PURGED ARRAYS OF TEMPERATURES 
AND DEPTHS AND FINDS THE DEPTHS OF ISOTHERMS SELECTED 
IN SUBROUTINE INTERP 

VARIABLES USED IN THIS SUBROUTINE 

TX,X— AS DEFINED IN SUBROUTINE INTERP 

IT1,IT2* IT3» IT4, IT5» IT6 — AS DEFINED IN SUBROUTINE 
INTERP 

Z11,Z10,Z09,Z08,Z07,Z06 — AS DEFINED IN MAIN PROGRAM 

KK— THE SUBSCRIPT OF THE DEPTH VALUE IN THE 2**N ARRAY 

DIMENSION TX<1),X< 1) ,Ziim,Z10(l ),Z09!1) ,Z08<1), 
*Z07<1),Z06<1) 

DO 1 1=1,100 

IF(TX(I).LE.IT1) GO TO 2 

1 CONTINUE 

2 Z 1 1 ( KK )=X( I) 

DO 3 J=I,100 

IF(TX( J) .LE.IT2) GO TO 4 

3 CONTINUE 

4 Z10I KK )=X( J) 

DO 5 K=J , 100 

IFITXIK) .LE.IT3) GO TO 6 

5 CONTINUE 

6 Z09( KK )=X IK) 

DO 7 L=K ,100 

IFITX(L) .LE.IT4) GO TO 8 

7 CONTINUE 

8 Z08(KK)=XIL) 

DO 9 M=L, 100 

IF(TXIM) .LE. IT5) GO TO 10 

9 CONTINUE 

10 Z07 I KK)=X ( M) 

DO 11 MM=M, 100 
IF(TX(MM)«LE*IT6) GO TO 12 

11 CONTINUE 

12 Z06 ( KK ) =X I MM ) 

RETURN 

END 



SUBROUTINE FINDZC ITX, ITC, KN, ICOUNT I 
THIS SUBROUTINE DETERMINES THE SUBSCRIPT OF THE DEPTH 
CORRESPONDING TO ITC AS FOUND IN SUBROUTINE INTERP 

VARIABLES USED IN THIS SUBROUTINE ARE 
TX,ITC — AS DEFINED IN SUBROUTINE INTERP 

ICOUNT — AS DEFINED IN SUBROUTINE PURGE 
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KN — THE SUBSCRIPT CF THE DEPTH ASSOCIATED WITH ITC 

D I MENS ION TX ( 1 ) 

DO I 1 = 1,1 COUNT 
I F ( TX ( I ) .LE. ITC) GO TO 11 
1 CONTINUE 
11 KN= I 
RETURN 
END 



SUBROUTINE NONGA { Z . T, I COUNT , ZT , ZC ,ZBT , KN ± JL ) 

THIS SUBROUTINE TAKES THE INPUT ARRAYS OF TEMPERATURE 
AND DEPTH, DIVIDES THE DEPTHS AT THE DEPTH CORRESPONDING 
TO ITC AND PRODUCES TWO SYMMETRICAL DEPTH ARRAYS, ONE A 
•REFLECTION* OF THE DEPTHS FROM Z(KN> TO THE SURFACE AND 
THE OTHER A 'REFLECTION* OF DEPTHS FROM Z(KN) TO THE 
BOTTOM CF THE TRACE. 

VARIABLES USED IN THIS SUBROUTINE ARE 
Z,T — AS DEFINEO IN SUBROUTINE INTERP 

ICOUNT — AS DEFINED IN SUBROUTINE PURGE 

ZT , ZC , ZBT — AS DEFINED IN MAIN PROGRAM 

JL — THE SUBSCRIPT OF THE DEPTHS ZT,ZC,ZBT IN THE ARRAY 
OF 2**N TERMS USED IN THE FOURIER ANALYSIS 

KN — THE SUBSCRIPT CORRESPONDING TO THE DEPTH OF THE 
ESTIMATED CENTER OF THE THERMOCLINE 

Z5 — THE 'REFLECTED* ARRAY OF DEPTHS FROM Z(KN) TO THE 
SURFACE 

Z3 — THE 'REFLECTED ARRAY OF DEPTHS FROM Z(KN) TO THE 
BOTTOM OF THE TRACE 

NUMBER tKDUMMY , NNN , NUMMER — COUNTERS USED IN THE 
PROCESS OF PRODUCING GAUSSIAN THERMOCLINES FROM THE 
NON-GAUSS I AN INPUT THROUGH 'REFLECTION* 

SK — A WEIGHTING FACTOR, WHICH IS APPLIED TO THE 
VARIANCE TO COMPENSATE FOR A LACK OF ACCURACY IN 
COMPUTATION DUE TO USE OF FINITE DIFFERENCE TECHNIQUES 
RATHER THAN INTEGRATION. ITS VALUE IS DEPENDENT UPON 
THE TEMPERATURE INTERVAL IDELTA T> USED IN COMPUTATION 

SUBROUTINES CALLEO FROM THIS SUBROUTINE ARE 
GAUSS1 
GAUSS2 

DIMENSION Z( 200) ,Z3( 200) ,Z5< 200 ) , T [ 200) , ZC C 1 ) ,ZT( 1) , 
♦ZBT ( 1 ) 

IF<(T< 1) — T 12) ) . GE . 10 . ) SK=l. 30000 
IFC(T{ 1 ) — T (2 ) ) . LE . 5. ) SK= 1.470000 
NUMB ER= I COUNT -KN + 1 
KDUMMY=ICOUNT-KN 
NNN=(2*{ ICOUNT-KN) ) + l 
ZCI JL)=Z(KN) 

DO 73 11=1, KN 
Z5 ( 1 1 ) = ZU1) 

73 CONTINUE 
I Z 1=KN- 1 
DO 7 12=1, IZ1 
X = Z ( KN ) - Z ( KN- 12) 

Z5IKN+I2) = Z ( KN ) + X 
7 CONTINUE 

NUMMER=2*KN— 1 
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CALL GAUSS1(Z5»NUMMER» I COUNT , SK,ZT,ZC ,JL) 
A KKK=KN-1 

DC 7 A 13=1, NUMBER 

Z3 ( KDUMMY + 13) = Z(KKK+I3) 

7 A CONTINUE 

DO 59 I A=l, KDUMMY 
X = Z ( KN+ 1 A ) - Z ( KN ) 

Z3 ( NUMBER - IA) = Z(KN) - X 
59 CONTINUE 

CALL GAUSS2(Z3,NNN, ICOUNT , SK, ZBT , JL) 

RETURN 

END 



SUBROUTINE GAUSS 1 <X,NT, ICOUNT , SK, ZT , ZC , JL ) 

THIS SUBROUTINE FINDS THE DEPTHS OF THE TOP AND CENTER 
OF THE THERMOCLINE THROUGH GAUSSIAN METHODS 

VARIABLES USED IN THIS SUBROUTINE ARE 

X — CORRESPONDS TO THE Z5 ARRAY OF SUBROUTINE NONGA 

NT— THE NUMBER OF ELEMENTS IN THE ARRAY OF X 

Zl— THE DEPTHS USEC IN COMPUTATION OF THE GAUSSIAN 
THERMOCLINE. BY DEFINITION, THE AVERAGE VALUE OF 
CONSECUTIVE DEPTHS IN THE X ARRAY, THE MIDPOINT 
OF 'DELTA T', AT WHICH THE DEPTHS MAKE THEIR 
CONTRIBUTION 

Z2 — THE ARRAY OF VALUES CORRESPONDING TO THE SQUARES 
OF THE VALUES OF Zl 

JL — AS DEFINED IN NONGA 

SK — AS DEFINED IN NONGA 

SIGSQR— THE VARIANCE CF THE ARRAY OF Zl 
SUM1 — THE SUM OF THE ARRAY OF Zl 

SUM2 — THE SUM OF THE ARRAY OF Z2, OR THE SUM OF THE 
SQUARES OF THE VALUES OF THE ARRAY OF Zl 

DIMENSION Zl(200),Z2(200),X(l),ZC(l),ZT(l) 

N=NT— 1 
SUM1 =0.0 
SUM2=0 .0 
DO 1 1=1, N 

zim«<xm*x(m>)/2. 

Z2 ( I )=Z1 ( I )**2 
SUM1=SUMH-ZI(I > 

1 SUM2=SUM2*Z2(I) 

ZC<JLJ=SUM1/N 

SI GSQR=( SUM2/NJ — ( ZC(JL)* t *2l 
ZT(JL)=ZC(JL )-SK*SQRT (SIGSQR) 

RETURN 

END 



SUBROUT INE GAUSS2 C X, NT .ICOUNT , SK, ZBT, JL ) 

THIS SUBROUTINE FINDS THE DEPTH OF THE BOTTOM OF THE 
THERMOCLINE THROUGH GAUSSIAN METHODS 

VARIABLES USED IN THIS SUBROUTINE ARE 

X — CORRESPONDS TO THE Z3 ARRAY OF SUBROUTINE NONGA 



67 



nooonnnnn nnnnnnnononnonn nnnnnnon 



NT — THE NUMBER OF ELEMENTS IN THE ARRAY OF X 

JlfSK — AS DEFINED IN SUBROUTINE NONGA 

SUMlfSUM2»Zl» 12 — AS DEFINED IN SUBROUTINE GAUSSI 

SIGSQR — AS DEFINED IN SUBROUTINE GAUSSI 

DIMENSION Z1 ( 200 ) »Z2( 200 ) » ZC ( 600 ) » X( 1),ZBT<1) 

N=NT— 1 

SUM1 =0.0 

SUM2=0 .0 

DO 1 I =1 , N 

Z 1 ( I )=<X( I )+X( 1 + 1) )/2. 

Z2(I )=Z1( 11**2 
SUM1=SUM1+Z1 ( I ) 

1 SUM2=SUM2+Z2 ( I I 
ZC I JL ) =SUM I/N 

SIGSQR=( SUM2/N)-(ZC( JL>**2) 

ZBT( JL )=ZC(JL )+SK*SQRT (SIGSQR ) 

RETURN 

END 



SUBROUTINE AVERGE ( AZ I N, AZOUT ) 

THIS SUBROUTINE TAKES THE ARRAY OF DEPTHS OF A SELECTED 
FEATURE (ISOTHERM OR THERMOCL I NE ) , AND PRODUCES AN ARRAY 
WITH A MEAN VALUE OF ZERO 

VARIABLES USED IN THIS SUBROUTINE ARE 

AZIN — THE ARRAY OF INPUT VARIABLES, WITH AN ARBITRARY 
MEAN VALUE 

AZOUT — THE OUTPUT ARRAY OF VARIABLES, WITH A MEAN OF 
ZERO 

ZINBAR— THE SUMMATION OF VALUES OF AZIN 

AZINB— THE MEAN OF VALUES OF AZIN 

DIMENSION AZIN(l) , AZOUT ( 1 ) 

ZINBAR=0.0 
DO 1 1=1,512 

1 ZINBAR=ZINBAR+AZIN( I ) 

AZINB=ZINBAR/512. 

DO 2 1=1,512 

2 AZOUT ( I )=AZIN ( I )- AZINB 
RETURN 

END 



SUBROUTINE HANN (ZIN,BOUT) 

THIS SUBROUTINE TAKES THE FOURIER COEFFICIENTS AS FOUND 
IN SUBROUTINE RHARM, AND HANNS THEM, TO REDUCE LEAKAGE 
BETWEEN NEIGHBORING FREQUENCY BANDS. 

VARIABLES USED IN THIS SUBROUTINE ARE 

ZIN — THE ARRAY OF FOURIER COEFFICIENTS FOUNO BY RHARM 

BOUT— THE HANNED FOURIER COEFFICIENTS 

DIMENSION ZI N ( 1 ) , BOUT ( 1 ) 

DO 1 1=1,2 

BOUT! I )= , 5*Z INC I )-.25*ZIN( 1 + 2) 

1 CONTINUE 

DO 2 1=509,510 
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BOUT! I+2)=-.25*ZlN( I >*.5*ZIN( 1+2) 

2 CONTINUE 

DO 3 1=3,510 

BOUT ( I )=- • 25*Z IN ( 1-2 )♦ • 5*Z INI I )-. 25*Z IN! I +2) 

3 CONTINUE 
RETURN 
END 



SUBROUTINE SMOOTH ( ARC I N, M, ITEST , ARGOUT , Z J , ICOUNT ) 

THIS SUBROUTINE TAKES THE HANNED FOURIER COEFFICIENTS 

AND AVERAGES THEM OVER FREQUENCY BANDS OF SPECIFIED 
WIDTH. THE NUMBER OF FREQUENCIES INCLUDED IN EACH BAND 
IS EQUAL TO M, AS READ IN BY THE MAIN PROGRAM BEFORE 
THIS SUBROUTINE IS CALLED. 

VARIABLES USED IN THIS SUBROUTINE ARE 

ARGIN — THE ARRAY OF HANNED FOURIER COEFFICIENTS 

M — A GROUPING FACTOR, AS DISCUSSED ABOVE 

ITEST — THE NUMBER OF GROUPS OF 2M FOURIER COEFFICIENTS 
WHICH IS ALSO ONE LESS THAN THE TOTAL NUMBER OF 
SMOOTHED FREQUENCY BANDS 

ARGOUT — THE ARRAY OF ITEST*! FREQUENCY BANDS 
ZJ — AS DEFINED IN MAIN PROGRAM 

ICGUNT — THE TOTAL NUMBER OF FREQUENCY BANDS IN ARGOUT 

FACTOR — A WEIGHTING FACTOR USED IN SMOOTHING 

N — THE NUMBER OF VALUES IN THE ARRAY OF DEPTHS 
SUBJECTED TO THE FOURIER ANALYSIS 

DIMENSION ARG iN(l>, ARGOUT! 1),ZJ( 1) 

CM=M 

N=512 

FACT OR = N/CM 
I TEST=N/{ 2*M) 

JTEST= I N— 1 )/M 

IF(ITEST*2.GE.JTEST) GO TO 3 
ITEST=ITEST+1 
3 CONTINUE 

IC0UNT=ITEST+1 
D01 KK=1» 512 
ARGOUT (KK) =0.0 
ARGIN! KKJ=ARGIN!KK)**2 

1 CONTINUE 
PART =0.0 

DO 2 J J=3 » M 

2 PART=PART+ARGIN! J J ) 

ARGOUT! 1)=FACT0R*( ARGIN! 1)*! 2. *PART> ) 

ZJ! l)=l. 

DO 100 11=2, ITEST 
MM=! 1 ! 2*II)-3)*M>+2 
MST0P=!!2*II)-1)*M+1 
DO 101 I J=MM, MSTOP 
101 ARGOUT ! 1 1 )= ARGOUT !II)+ARGIN!IJ) 

ARGOUT! II )=F AC TOR* ARGOUT (II) 

130 ZJ(in=II 
MN=MM+M 

ZJ! ICOUNT )=ICOUNT 
IDIFF=5Ii-(( ITEST*2>-1>*M 
IF! 2*! IDIFF/2) .EQ.IDIFF) GO TO 30 
IDIFF=IDIFF*1 
30 PRDIFF= IDIFF/2 
DIFF=512./PR0IFF 
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00 102 JK=MN *512 

102 A R GOUT ( I COUNT ) = ARG0UT ( I COUNT ) +ARGINIJK) 
ARGOUT ( I COUNT ) =DI FF*ARG0UT ( I COUNT ) 
RETURN 
END 



SUBROUTINE PER IOD ( M , I TEST , APER , A J , B J ) 

THIS SUBROUTINE IS USED TO GENERATE THE PERIODS 
ASSOCIATED WITH THE SMOOTHED SPECTRAL BANDWIDTHS, AND TO 
GENERATE THE FINAL SET OF COUNTERS FOR USE IN GRAPHICAL 
DISPLAY 

VARIABLES USED IN THIS SUBROUTINE ARE 
M — AS DEFINED IN SUBROUTINE SMOOTH 

ITEST — AS DEFINED IN SUBROUTINE SMOOTH 

A PER — THE ARRAY OF PERIODS ASSOCIATED WITH THE ARRAY 
OF SMOOTHEC BANDWIDTHS 



AJ — AS DEFINED IN MAIN PROGRAM 

BJ — DEF INED AS THE NATURAL LOGARITHM OF AJ FOR USE IN 
PRODUCING A SEMILOG GRAPHICAL OUTPUT 

DIMENSION APER(l) » AJ( 1) »BJ( 1 ) 

I PTS= I TEST + l 
APER(1)=(2.*512. )/M 
DO 1 1=2 1 1 TE ST 
APERt I >=512. /( ( 1-1 )*M) 

APERt IPTS)=2048./t 1024.- ( 512. -< t 2* I TEST+1 ) *M ) ) ) 

DO 3 1 = 1 » I PT S 

AJCI )=512./<M*APER(I ) ) 

DO 2 1=1, I PT S 
B J ( I ) =ALOG ( A J ( I ) ) 

RETURN 

END 



SUBROUTINE GROUP t A J ,BJ , APER, I PTS , ENDT , ENDC , ENDB, END1 1 , 
♦END10, END09, END08 , END07, END06 > 

THIS SUBROUTINE IS USED TC AVERAGE THE RAW SPECTRA 
PRODUCED IN 'IGROUP' TIMES THROUGH THE MAIN PROGRAM , 

AND TO PRINT THE VALUES ASSOCIATED WITH THESE AVERAGED 
SPECTRA 

VARIABLES USED IN THIS SUBROUTINE ARE 

A J ,BJ , APER — AS DEFINED IN SUBROUTINE PERIOD 

I PTS — THE TOTAL NUMBER OF SPECTRAL BANDWIOTHS COMPUTED 
IN SUBROUTINE SMOOTH 

END11,END10, ENDT, ENDC, ETC. --AS DEFINED IN MAIN PROGRAM 

DIMENSION AJ( 1 >,BJ( 1) , APERt 1> , ENDT Cl) , ENDC ( 1 ), ENDB (1) , 
*END1 1(1), END 10(1), END09(1),END08(1), END07 ( 1) ,END06( 1) 
DO 1 1 = 1, I PT S 

ENOTt I )= ,5*( ENDT ( II + ENDTI I+IPTS)) 

ENDC ( I )=. 5*1 ENDC < I I+ENDCI I+IPTS) ) 

ENDBI I )=.5*t ENDB( I )+ENDB( I+IPTS)) 

EN01K I)=.5*(END11( I ) +END1 1 ( I+IPTS) ) 

ENDlOt I )=.5*t ENDlOt I )+END10t I+IPTS)) 

EN009C I ) = • 5*( END09I I ) +END09I I+IPTS)) 

END0B(I)=.5*( ENDO 8 ( I)+END08( I + IPTS)) 

END07 ( I ) = .5*t ENDO 7(1) +END071 I+IPTS) ) 
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END06 ( I )= . 5* ( END06 ( I ) *END06( I + IPTS) ) 

1 CONTINUE 
WRITE! 6,30) 

30 FORMAT 1 1H1 ) 

WR I TE ( 6 , 3 ) 

3 FORMAT ( / /40X » 'AVERAGED SPECTRAL ESTIMATES',/) 
WRITE(6,4) 

4 FORMAT I3CX,' NR * ,4X,' PERIOD' ,IOX,'TOP' ,4X, 'CENTER' ,4X, 
♦'BOTTOM',/) 

DO 5 I =1 , I PTS 

5 WRITE (6,6) AJI I) , APERl I ) ,ENOT (I ) ,ENDCC I ) ,ENOBI I ) 

6 FORMAT (28X.F4. 1,4X, F6 . 2 , 4X , 3F 10 . 2 ) 

WR I TE ( 6,30) 

WRITE (6,3) 

WR ITE(6»7) 

7 F0RMAT(30X, 'NR' ,4X, 'PERIOD* ,7X, * II DEG*,4X,'10 DEG* , 

♦ 5X » '9 DEG' ./) 

DO 8 I = 1 , I PT S 

8 WRITE (6,6) AJ( I) , APER ( I ) , END 11(1) ,ENDI0( I ) » END09 ( I) 

WR IT E ( 6, 30 ) 

WRITE(6,3) 

WRITF ( 6, 10) 

10 FORMAT OCX, 'NR', 4X, 'PERIOD* ,8X, *8 DEG',5X,'7 DEG',5X, 
♦'6 DEG' ,/) 

DO 9 1 = 1 , I PT S 

9 WRIT E (6,6 ) AJ( I ),APER( I ),END08( I ),END07( I )»£ND06( I) 
RETURN 

END 
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